Mon. Not. R. Astron. Soc. 000,[T][T7](2012) Printed 15 November 2012 (MN KOgX style file v2.2) 



Angular momentum transfer to a Milky Way disk at high redshift 

H. Tillson 1 *, J. Devriendt 1 ' 2 , A. Slyz 1 , L. Miller 1 & C. Pichon 1 ' 3 

1 Department of Physics, University of Oxford, The Denys Wilkinson Building, Keble Road, Oxford, 0X1 3RH, UK 

2 Centre de Recherche Astrophysique de Lyon, UMR 5574, 9 Avenue Charles Andre, F69561 Saint Genis Laval, France 
s Institut d' Astrophysique de Paris, 98 his boulevard Arago, F-75014 Paris Cedex, France 



15 November 2012 

ABSTRACT 

An Adaptive Mesh Refinement cosmological resimulation is analyzed in order to test whether 
filamentary flows of cold gas are responsible for the build-up of angular momentum within 
a Milky Way like disk at z ^ 3. A set of algorithms is presented that takes advantage of the 
high spatial resolution of the simulation (12 pc) to identify: (i) the central gas disk and its 
plane of orientation; (ii) the complex individual filament trajectories that connect to the disk, 
and; (iii) the infalling satellites. The results show that two filaments at z > 5.5, which later 
merge to form a single filament at z < 4, drive the angular momentum and mass budget of the 
disk throughout its evolution, whereas luminous satellite mergers make negligible fractional 
contributions. Combined with the ubiquitous presence of such filaments in all large-scale cos- 
mological simulations that include hydrodynamics, these findings provide strong quantitative 
evidence that the growth of thin disks in haloes with masses below 10 12 Mq, which host the 
vast majority of galaxies, is supported via inflowing streams of cold gas at intermediate and 
high redshifts. 

Key words: galaxies revolution - galaxies: formation - galaxies:haloes - galaxies:high red- 
shift - methods:numerical 



1 INTRODUCTION 

Our understanding of the origin of angular momentum within 
galaxies dates back to the pioneering works of Hoyle ( 195TJ and 
Peebles] (1969[ >, who modelled protogalaxies as spherical Eule- 
rian patches and argued that torques are exerted on a given patch 
due to tidal interactions with neighbouring patches. Doros hkevich| 
( 1970 1 removed the assumption of spherical symmetry and demon- 
strated that as perturbations follow the expansion of the back- 
ground, their total angular momentum evolves with the scale factor 
as J [a) oc a 3 ' 2 in an Einstein-de Sitter universe. An iV-body sim- 
ulation study by White ( 1984) later confirmed that this growth of 
angular momentum continues until the perturbation, modelled in a 
Lagrangian framework, breaks away from the background expan- 
sion and collapses to form a halo, whereupon its angular momen- 
tum remains constant. White ( 1984) hence argued that the magni- 
tude of the angular momentum acquired for a given halo is set at 
the epoch of maximum expansion. 

A decade after the seminal papers by |Peebles| ( |1969| > and 
|Doroshkevichl p970) , pall & Efstathiou| J 1980) expanded beyond 
the early predictions of tidal torque theory to account for the hier- 
archical growth of dark matter haloes described in White & Rees| 
( 1978 ), and consequently painted a simple picture of galaxy disk 
formation. These authors claimed that gas residing in the poten- 
tial wells of collapsed density perturbations is shock heated to the 
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virial temperature of the halo, and that the inner gas regions subse- 
quently cool and lose their pressure, sinking to the halo centre on a 
free-fall timescale. Galaxy disks hence form from the 'inside-out', 
as the infalling gas retains its specific angular momentum. Several 
authors have since reported that the mean stellar age decreases ( jde] 
|Jong|1996[ |MacArthur et al.|2004| |Gogarten et al.|2010|> and s pe- 



cific star formation rates increase ( Munoz-Mateos et al. 2007 ) as 
a function of radial position from the disk centre, thereby lending 
support to the inside-out paradigm of disk growth. 

This standard model of disk galaxy formation makes two im- 
portant assumptions: (i) the specific angular momentum distribu- 
tions of the gas and virialised dark matter halo are initially equal, 
and; (ii) the specific angular momentum of the gas is conserved 
upon collapse to the central disk. This simple model has two ap- 
pealing features. Firstly, by making a prediction for the relationship 
between the disk scale length and the host halo virial radius, it is 
able to recover locally observed distributions of disk scale lengths 
and the J -band Tully-Fisher relation (e.g. |Dalcanton et al.|1997| 
|Mo et al.|1998] |de Jong & Lacey|2000| l. The second success for 
this model lies in the discovery, using numerical simulations, of the 
remarkable property that the cumulative fraction of mass within a 
halo that has specific angular momentum less than j, denoted by 
M(< j), is well fit by a universal function that depends on a sin- 
gle free parameter (Bullock et al. 2001 I. This result holds for both 
dark matter and (non-radiative) gas present in haloes of all masses, 
although the angular momentum vectors of these two components 
are not generally aligned (van den Bosch et al. 2002). 
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Despite these successes, recent years have witnessed a resur- 
gence in tackling one of the defining characteristics of the standard 
theory of disk formation: that gas crossing a halo's virial sphere is 
shock heated to the virial temperature of the halo, which typically 



corresponds to X-ray temperatures of a few 10 6 K. Binney | jl977f 



was the first to claim that, contrary to the standard paradigm, high 
density gas in low mass haloes (i.e. most haloes at high redshift) 
probably crosses an isothermal shock front close to the disk rather 
than an adiabatic one at the host halo virial radius, and hence ar- 
gued that the vast majority of gas is expected to stream in towards 
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the disk in cold flows with temperatures close to 10 K 
& Dekel (2003) reignited this idea by demonstrating that for haloes 
less massive than ~10 11,6 Mq, gas cools efficiently, loses some of 
its pressure, and is unable to support a virial shock. It was hence 
concluded that for haloes of mass M H < 10 116 Mq, gas flows 
to the disk in a cold mode as opposed to a hot mode. The lack 
of redshift evolution in this threshold mass and the penetration of 
cold streams deep inside the virial region of most haloes was later 
confirmed by several numerical studies ( Keres et al. 2005 , Ocvirk| 
|et al.|20 08). The existence of a bimodal gas accretion phase poten- 
tially has enormous implications for star formation, and may. as the 
above authors have speculated, account for the observed lack of soft 
X-ray flux and copious Lyman-a emission in high redshift galax- 
ies l |Birnboim & D ekel 2003 ), the bimodality in galaxy colours and 
its relation to galaxy morphology, including the existence of red 
galaxies at z > 1 and massive bursts of star formation at z ~ 2-4 
( |Dekel & Birnboim|2006), and possibly star formation rate down- 
sizing 1 Kere s et al.|2005||Ocvirk et al.|2008) . These are hence very 
exciting times for galaxy disk theory. 

Having established that gas inflow in most haloes at high red- 
shift proceeds via cold flows that often terminate in the outskirts of 
gaseous disks ( |Brooks et al.|2009[ |Powell et al.|2011| >, one natu- 
rally wonders about the processes governing the transport of angu- 
lar momentum onto disks within these overdense regions. This is 
undoubtedly a complex, multiscale problem. On the largest scales 
of several Gpc, the Universe is thought to be arranged in a complex 
web of cosmic structure ( |Bond et al.| 1 996 ; Pogo syan et al.|1998| >, 
and simulations with differing spatial resolution have helped to de- 
compose the intricate patterns of this web (e.g. Sousbie 20TTJ. It 
appears that large sheets surrounding Gpc-scale voids intersect to 
form filaments that funnel gas from the voids to the intersection 
nodes of the web, where haloes of dark matter form. The amount 
of mass shared between these separate phases has recently been 
the subject of debate and involves examining the tidal field tensor 
that describes the second order derivatives of the gravitational field, 
in bot h the linear ([ DorosJikevichM^70 1 and non-linear ( |Shen et aT] 
|2006||Hahn et al.|2007||Arag6n-Calvo et al.|2010| l regimes of per- 
turbation growth. Yet the general conclusions of these studies are in 
agreement: filaments dominate the mass budget. It is hence likely 
that this component also carries large amounts of angular momen- 
tum to the halo nodes of the cosmic web. 

Several recent studies have shed some light on the sub- 
ject of angular momentum transport by filaments on galactic 
scales. Da novich et al.|p012| > performed a statistical study of 350 
Milky Way size dark matter haloes selected from the HORIZON- 
MARENOSTRUM simulation (details given in |Ocvirk et aT1|2008| 
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mass Mh — 10 Mq, whose luminous components were sim- 
ulated with a physical resolution of lkpc. They found that the 
streams of cold gas flowing towards the disk were oriented in a 
narrow plane, and that the angular momentum transported by this 
infalling component was highly misaligned with respect to the an- 



gular momentum direction of the disk, until the approximate disk 
boundary, whereupon it dramatically swung into close alignment. 
An earlier study by |Pichon et al.| ( f2011| > analyzed several outputs 
from the same run and examined the nature of filament trajecto- 
ries on kpc scales, in an attempt to understand the existence of thin 
gas disks at high redshift that are thought to have formed from the 
inside-out. These authors demonstrated that material accreted at the 
virial sphere carried more angular momentum at later times, and 
attributed this phenomenon to a 'lever' mechanism: recently ac- 
creted gas has a larger impact parameter owing to the velocity sway 
of filaments on large scales. It was further hypothesized that these 
large-scale drift velocities arise from the asymmetric cancellation 
of motions of gas pumped out of voids. [Fiction et af] ((201 1 ) hence 
concluded that the angular momentum transported along cold gas 
flows into halo virial regions originates from the asymmetry of sur- 
rounding large-scale voids. Kimm et al. (20111 investigated these 
claims by probing an individual halo to a higher physical resolu- 
tion of ~50 pc, and found that within the virial sphere, the specific 
angular momentum of the radiative gas was systematically larger 
than that of the dark matter. These authors argued that this was a 
manifestation of a dark matter angular momentum diffusion caused 
by the mixing of dark matter particles carrying different amounts 
of angular momentum at the time of accretion (this mixing is due to 
the dark matter undergoing shell-crossing as it passes through walls 
and filaments, before entering the host's virial region). They justi- 
fied this claim by demonstrating that there was a large distribution 
in the age of the dark matter particles at any given radius (see also 
the impact parameter distributions and backslash velocities mea- 
sured at the virial radius in Figs. 15 and 20 of Aub ert & Pichon| 
|2007[ >. On the other hand, the halo gas does not experience the 
same angular momentum diffusion since it streams almost exclu- 
sively along filaments (i.e. in a preferential direction) towards the 
central region where it eventually experiences an isothermal shock 
and remains trappe d jBirnboim~& Dekel 2003; Keres ^Tal.|2005| 
|Ocvirk etal . 2008 1. Kimm et al.| ( |201 l| l demonstrated that this does 
indeed happen by showing that the amount of specific angular mo- 
mentum transported by gas was constant with radius at both low 
(0 ^ z ^ 3) and high (z > 3) redshifts, except for r < 0.1r v i r , 
whereupon it fell dramatically. This hints at unresolved complex 
dynamics within the 'disk' region (e.g. |Book et al.|201 1| | of high 
redshift galaxies. 

Yet due to physical resolution constraints, the |Kimm et al.| 
l |201 l) , |Pichon et al.| l [20TT) and |Danovich et al.| < |20l2l > studies were 
neither able to accurately resolve the disk scaleheight and scale- 
length at high redshift, nor the dynamics within the disk region. By 
analyzing outputs from one of the high resolution NUT simulations 
(see |Powell et al.| "201 1 I, which tracks the high redshift evolution of 
a single Milky Way like galaxy with a physical resolution scale of 
12 pc, this paper rises to the challenge of attempting to account for 
the amount of angular momentum locked-up in a Milky Way like 
disk that resides in a halo fed by streams of cold gas, and comple- 
ments the papers listed above. The aim is to ascertain whether fila- 
ments dominate the angular momentum budget of the central disk 
across time. This would provide quantitative evidence in favour of 
the cold mode of disk growth, a theory currently lacking obser- 
vational con firmation JFaucher-Giguere et al.|2010| [Steidel et al.| 
[20T0l but see |Fumagalli et al.|2011) ~ 

This paper is arranged as follows. In Section [2] the details 
of the simulation and the physics encoded within it are briefly 
discussed. The various algorithms that have been developed to: 
a) compute the orientation plane of the resimulated Milky Way 
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like disk; b) resolve the individual gas filaments and satellites that 
merge onto it, and; c) quantify the angular momentum locked-up in 
each of these galaxy components, are then introduced in Sections 
[3] and [4] A more detailed description of the filament and satellite 
identification algorithms is provided in Appendices [A] and [B] re- 
spectively. Section[3]presents the results, Section[6]discusses their 
implications, and Section[7]summarizes our conclusions. 



2 THE SIMULATION 
2.1 The NUT cooling run 

We examine one of the NUT simulations to understand how a Milky 
Way like disk acquires its angular momentum. The NUT suite (Pow- 
|ell et al.|2011) is a set of ultra-high resolution simulations that em- 
ploy a zoom technique (see e.g. |Navarro et al.| 1995^ to resimulate 
a Milky Way like galaxy and follow its evolution across redshift in 
a ACDM cosmology. The resolution and number of physical pro- 
cesses that are modelled vary from one simulation to the next. The 
highest resolution runs in the NUT suite have a maximum spatial 
resolution Ai res of 0.5 pc but terminate at high redshifts (z ~ 6) 
due to computational time restrictions, and so in this paper a lower 
resolution Aa; res = 12 pc run is analyzed as it probes redshifts 
down to z = 3. More specifically, the NUT cooling run, hereafter 
NutCO run, examined in this paper has been performed using the 
publicly available Adaptive Mesh Refinement (AMR) code RAM- 
SES {Tteyssier 20021, incorporating the following physics: 

(i) Epoch of reionization — a spatially uniform, redshift- 
dependent ultraviolet (UV) background instantaneously switched 
on at z = 8.5 (following the |Haardt & Madau|1996| UV model). 

(ii) Cooling and star formation. 

The important technical details of the NutCO run are now 
briefly summarized. The chosen resimulated Milky Way like halo 
satisfies two conditions: it is a typical 2a density fluctuation, i.e. it 
is not embedded in a denser group or cluster environment at z — 0, 
and its present day mass (Mh — 5 x 10 11 Mq) is com parable 
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to the threshold mass (M H ~4x 10 11 Mq) reported by 
et al. (20081 below which gas streams along filamentary structures 
at low temperatures of T < 2 x 10 4 K. The entire simulation box has 
a comoving side length of 9 h^ 1 Mpc and starts at z = 499, with 
the initial conditions generated using the package MPgrafic ( |Prunet| 
|et al.|20 08). The resimulation box around the main halo has a co- 
moving side length of ~2.7ft _1 Mpc, and the simulation evolves 
in a universe of WMAP5 cosmology (Dunkley et al. 2009| with 
n M = 0.258, n A = 0.742, Q b = 0.045, a s = 0.8 and h = 0.72. 
The coarse root grid for the entire simulation has 128 cells along 
each dimension of the 9/i _1 Mpc box, whereas the resimulated 
region contains three higher resolution nested grids yielding an 
equivalent particle resolution of 1024 3 dark matter particles, each 
with mass Mum = 5.4 x 10 4 Mq. A quasi-Lagrangian refine- 
ment scheme is used to maintain Ax rcs = 12 pc in physical co- 
ordinates as the simulation evolves. Higher levels of refinement 
are spawned once either: a) the baryonic mass in the cell exceeds 
8?tisph (mspH = 9.4 x 10 3 Mq) or; b) the number of dark matter 
particles in the cell exceeds eight. This scheme therefore strives to 
achieve a roughly equal gas mass per cell. 

The NUT simulation includes stars and so a brief discus- 
sion of the prescriptions used for modelling star formation, which 
are described in detail by Rasera & Teyssier (20061, is now pro- 
vided. The gas in a cell can cool to temperatures T ~ 10 4 K via 



bremsstrahlung radiation (effective until T ~ 10 6 K), and via col- 
lisional and ionization excitation followed by recombination (dom- 
inant for 10 4 < T/K sC 10 6 ). The metallicity of the gas is fixed to 
10~ 3 Zq to allow further cooling to molecular cloud temperatures 
of ~1K via metal line emission, and once the gas is sufficiently 
dense (i.e. p > po, where po is a density threshold and is chosen 
to represent the density of the interstellar medium, hereafter ISM), 
star formation naturally ensues, a process which is modelled by a 
Schmidt law (Schmidt 19591. The star formation efficiency param- 
eter e is fixed at 1% per free-fall time in concordance with observa- 
tions pCrumholz & Tan|2007| l and po is chosen to be equal to 400 
hydrogen atoms per cubic centimetre. If the Jean's scale Aj is close 
to or below the minimum resolved scale Ai m i„ (determined by the 
highest level of refinement ^ ma x), it is possible that gas clouds ar- 
tificially fragment, leading to artificial star formation. Therefore, 
when in excess of the ISM density po, the gas is forced to follow 
a polytropic equation of state. The number N of stellar particles 
that form in a cell whose gas is sufficiently cool and dense is drawn 
from a Poisson distribution, based on the cell mass, the local star 
formation timescale, and the minimum stellar mass allowed for a 
star particle m*, m i n , which itself is proportional to the minimum 
cell volume. In order to prevent an overproduction of stellar parti- 
cles in a given cell, all of the N newly formed stellar particles are 
then merged into a single particle of mass m* = 7Vm* >m i n . Stel- 
lar particles in the NutCO run analyzed in this study are therefore 
stellar clumps of up to several thousand stars. 



2.2 Identifying (sub)haloes & galaxies and following their 
merger histories 

The Most massive Sub-node Method subhalo-finding algorithm, 
hereafter MSM (Tweed et al. 2009), was used in order to detect 
the haloes, subhaloes, stellar clumps and stellar subclumps at each 
time output of the NutCO run. MSM assigns a local density esti- 
mate to each particle computed using the standard Smoothed Par- 
ticle Hydrodynamics kernel (Monaghan & Lattanzio 19851, which 
weights the mass contributions from the JV closest neighbouring 
particles (N — 20, 32 for the haloes and stellar clumps, respec- 
tively). Haloes are then resolved by imposing a density threshold 
criterion and by measuring local density gradients. In order to de- 
tect host haloes, stellar clumps and associated levels of substruc- 
ture, MSM successively raises the density thresholds on the host un- 
til all of its node structure has been resolved. The most massive leaf 
is then collapsed along the node tree structure to resolve the main 
halo or the main stellar clump and the same process is repeated for 
the lower mass leaves, defining the embedded substructures. 

Each resolved halo (stellar clump) object is forced to contain 
at least 40 (100) particles to ensure reliable detections. Note that 
the stellar clump resolution limit is set above the dark halo limit 
to avoid identifying groups of stars within the disk as independent 
objects. The TreeMaker code jTweed et al. 2009) is then used to 
link together all the time outputs by finding the fathers and sons of 
every halo and subhalo, and every stellar clump and subclump. 



3 RESOLVING THE GALAXY COMPONENTS 

This section describes the algorithms that have been used to iden- 
tify the likely sources of the disk's angular momentum within the 
host's virial region. The following terms — r v i r , disk, satellites and 
virial sphere — are henceforth made in reference to the NUT host 
halo. 
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Figure 1. Images of the filament gas within a box centred on the main host at 2 ~ 10 (left) and z ~ 8 (right) in the NutCO ran, gener- 
ated using the Visualization and Analysis Platform for Ocean, Atmosphere, and Solar Researchers software (commonly referred to as VAPOR, see 
https : / /www . vapor .ucar.edu/ for details). The sidelength of the box in the left panel corresponds to 4r vlI (~17 kpc) and this face-on view il- 
lustrates the clear separation of the filaments into three distinct regions (labelled accordingly) at large scales and early times. The image in the right panel 
shows the inner virial region (~7.5 kpc) at lower redshift when the dynamics of the system have changed, and the velocity flow tubes demonstrate the complex 
trajectories of gas flow within the gas disk region, colour coded according to the filament that they are likely to trace. The density of the gas has been dimmed 
and the viewing angle rotated to best highlight each trajectory. Both panels have been refined to £ = 16, corresponding to a physical cell width of ~17pc 
(left) and ~22 pc (right), and the colour bar indicates the gas density in dimensionless units of logio(n/cm — 3 ), where n is the number density of hydrogen 
atoms per cubic centimetre. Filament gas has been detected by applying the nT criteria, but gas with density logio(n/cm — 3 ) < — 1 has been excluded from 
both figures for clarity. The rightward red, upward green and forward blue axes correspond to the x, y and z directions, and the vertical and horizontal black 
planes in the left panel show the spatial cuts in the x and y directions that have been used to start the colouring of the tracer particles, which is a technique that 
has been developed to separate each of the filaments at high redshift (see Section [33] and AppendixjAj. 



3.1 Grid notations 

This study analyzes cubic and spherical grids Q x centred on the 
galaxy (i.e. on its most dense cell) with physical half-lengths and 
radii equal to a;r v i r , where x £ {0.1, 0.15, 0.5, 1.0, 2.0}. The nu- 
merical subscript x will be explicitly stated when discussing a cer- 
tain grid type, but for comparisons between grid types (e.g. the gas 
grid Q g or the stellar grid <?*) it is replaced by a symbol denoting 
the type, and the extent of the grid is implicitly understood. All 
of the spherical (cubic) grids in this paper with physical diameter 
(full-length) L are subdivided into n cubic cells of equal physical 
sidelength Ax and are hence 'fixed', obeying the relation: 



Ax 



(1) 



Since there are 2 cells along each dimension of the physical sim- 
ulation grid length L a (z) at level £, it follows that the refinement 
level of Q x is given by: 



\n{L s /Ax) 
In 2 



(2) 



where L s (z) — 9 h^ 1 Mpc/ (1 + z) and i is rounded to its nearest 
integer value. The routines that generate fixed density, temperature 
and velocity grids at a given location by performing an interpola- 
tion from the simulation grid for gas properties and a Cloud-in-Cell 
smoothing technique for particle properties, require L and £ as in- 
puts at runtime. Memory and CPU time constraints mean that n is 
not allowed to exceed 1024, however, so we decrease £ successively 
in integer units in equation |2]( until this constraint is satisfied. As a 
result, the resolution Ax for the grids enclosing the r v i r and 2r v i r 
regions deviates from the maximum resolution Ax TCS for low red- 



shift outputs, as the virial radius of the host (and hence L) expands. 
However, this 'time bias' does not affect the grids (Jo.i and Q0.15 
used for performing the angular momentum computations (see Sec- 
tion|4j, as their cell width remains equal to Ax res at all times. 



3.2 Disentangling the filaments using a tracer particle 
colouring algorithm 



|Powell et al.H20lT) have shown that parcels of gas whose number 
density of hydrogen atoms n and temperature T simultaneously 
satisfy 0.1 sC n/cm~ 3 ^ 10 and T/K sC 2 x 10 4 , are repre- 
sentative of gas belonging to filaments in the NUT simulations at 
high redshift (z ^ 9). The analysis in this work resolves the NUT 
filaments down to z = 3, hence in order to tackle the possible 
time evolution in their density, this paper uses the prescription be- 
tween the lower threshold and the background density introduced 
by |Kimm et al.|p01 1| >, who argued that large-scale filament density 
is gravitationally coupled to the expansion of the Universe. They 
defined the lower number density bound ul as: 



riL = 



SfpfbXu 



m H 



(3) 



where p, Sf, fb,Xn and run are respectively the mean density of 
the background, the density contrast of the filament (measured with 
respect to p), the fraction of the total mass in baryons (= Qt/^M), 
the primordial relative mass abundance of hydrogen (76%), and the 
mass of a hydrogen atom. A typical value of wl from the NutCO 
run at z = 3 is ~0.01 cm -3 . 

It was also found that an upper temperature limit of T u ^ 
2 x 10 5 K (e.g. |Keres et al.||2009| [Faucher-Giguere et al.||2011 l 
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Figure 2. These gas density images are z projections of the large-scale region surrounding the central disk in the NutCO run, and illustrate the spatial 
distribution of particles tracing the filaments at z ~ 9 (left panel) and 2 = 3 (right panel). The inner and outer solid circles have radii equal to r v - lr and 2r v - lT 
respectively, where r v i r ~ 5.3 kpc (left) and r v ; r ~ 31 kpc (right). The dashed circles represent the 2r v i r regions at the previous time output, corresponding 
to z ~ 9.5 (rvi r ~ 4.7 kpc) and z ~ 3.04 (r v i r ~ 30 kpc). All of the circles are centred at the location of the densest cell in the disk region at the relevant time 
output. Left panel: Those tracer particles that have been identified with each of the purple, yellow and red filaments at z ~ 9.5 and that coincide with filament 
cells refined to I = 15 (Ax ~ 38 pc) at z ~ 9, have been updated to their new positions. The offset in centres between the outer solid and dashed circles 
arises due to the motion of the centre of disk rotation, and the gap between each set of coloured tracers and the dashed circle highlights the displacement of 
the tracers over the time interval. Uncoloured new tracer particles that have crossed the 2r v ; r sphere over the interval and that move along filament trajectories 
at 2 ~ 9 have been omitted. Right panel: The filament configuration at i = 14 (Ax ~ 190 pc) for the final time output at z = 3. This image includes both 
updated coloured tracers and new tracers that have been assigned a colour using the nearest neighbours method (see text). In Section[3]it is demonstrated that 
the purple and yellow filaments start to merge around z ~ 5.5, with the former filament surviving, and this is confirmed by the absence of yellow tracers in 
this panel. Only 5% of all the coloured filament tracers are shown at this time output, for clarity. 




detected filaments to a higher level of accuracy for z ^ 3 in the 
NutCO run than the classical cold mode temperature of T u ^ 2 x 

10 4 



K (e.g. |Kay et al.|200o{ |Birnboim & Dekel|2003[ |Keres et al. 



2005 1 at which gas is thought to cool via Lyman-a emission (e.g. 
Fardal et al. 20011. There are two reasons for this. Firstly, once 
reionization has completed (z rc = 8.5), the UV background heats 
all the gas uniformly, without taking into account its capacity to 
self-shield from UV photons in denser environments. Secondly, the 
absence of metal production sources in the NutCO run causes the 
gas metallicity to remain fixed at its very low initial value (10 -3 
Zq). As a consequence, cooling is inefficient and filaments stay 
hotter for longer (e.g. Sutherland & Dopita 1993| >. 

Hence, in order to identify the neutral hydrogen and helium 
gas in the filamentary phase in the NutCO run, the density and 
temperature criteria, hereafter nT, reported by Powell et al.|p0TT) 
have been used for z 3, but the lower density threshold has been 
replaced by the value given by equation {3} and the upper tempera- 
ture limit has been increased to T u — 2 x 10 5 K. The left panel in 
Fig.[T]shows an image of the three filaments within 2r v i r of the host 
centre at z ~ 10 in the NutCO run, found by imposing these nT 
criteria. The filaments occupy three distinct regions on large scales 
at this epoch and the flow appears to be ordered and mostly radial 
(a discussion of how gas that is shared by filaments and satellites 
is apportioned between these two phases is preserved for Section 
3.3[ >. The dynamics of the gaseous motions in the central region 
are far more complex, however, as demonstrated in the right panel, 
which shows the velocity flow vectors of filament gas within the 
virial region at a later time (corresponding to z ~ 8) when the 
system has evolved to a different configuration. The individual fila- 
ment trajectories bend, twist and mix with each other in the central 
region, yet often remain intact. In order to test the hypothesis that 



the disk acquires its angular momentum from filamentary motions, 
it is critical to follow these separate individual filament trajectories. 
A suitable tool for performing this Lagrangian exercise is the use of 
tracer particles (e.g. Dubois et al. 20121, which have been demon- 
strated to preferentially trace filament flows in galaxy simulations 
(Picho n et al.|201l) . The tracer particles in the NUT suite have an 
identical spatial distribution to dark matter at the beginning of the 
simulation, and are assigned zero mass. They are simply advected 
with the gas and can therefore be used to track its velocity field 
within the host virial region. 

The algorithm that resolves each filament at high redshift is 
henceforth referred to as the 'tracer propagation method', and is 
described in detail in Appendix [A] which also includes a discus- 
sion of the (negligible) uncertainties associated with the method. 
The main idea is to assign the tracer particles within each fila- 
ment a unique colour at some early time, and then use a nearest 
neighbours scheme to construct the individual filament trajectories 
within the host virial region at later times. Throughout this paper 
reference is made to the purple, yellow and red filaments in regions 
1, 2, and 3 of the left panel of Fig.[T]respectively. The left panel of 
Fig.|2]shows tracer particles colour coded according to the filament 
they sample at z ~ 9.5, updated to their new filament positions at 
z ~ 9. The displaced tracers extend across the full virial region at 
z ~ 9 (given by the inner solid circle), and so new tracers enter- 
ing the 2r v h- region at this epoch (denoted by the outer solid circle) 
are appropriately coloured. Since colour is propagated across time, 
one is able to test the overall performance of the method by con- 
sidering the final output. Coloured filament tracers within 2r v i r at 
z = 3 are hence shown in the right panel of Fig. [5] (only 5% of 
all the coloured filament tracers are plotted for clarity, but this sub- 
sample is representative of the whole distribution). The purple and 
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Figure 3. Projections along the z direction of all the gas within a box of sidelength equal to 2r v - lr , corresponding to the diameter of the orange circles. The 
centre of each image coincides with the centre of the stellar clump of the host galaxy at z ~ 8 (left) and z = 3 (right), with i = 16 (Ace ~ 22 pc) and 
I = 15 (Ax ~ 95 pc). The virial regions of the stellar and dark halo satellites of the main host, as identified by the satellite-finding algorithm, are shown by 
the red and white circles respectively, with each circle centred at the densest member star identified as part of the satellite. The constituent satellite stars have 
been detected by the clump-finder, and are shown as small filled black circles. The white horizontal bars and vertical colour bars indicate the length scale (in 
physical kpc) and the density colour scale (in dimensionless units of logio number of hydrogen atoms per cubic centimetre) respectively. Note that the large 
stellar satellites in the upper right corners of both panels have not been identified as luminous satellites because their stellar virial regions do not infringe the 
host's virial sphere at these epochs. 



yellow filaments undergo a merger around z ~ 5.5 in the simula- 
tion (this is discussed in more detail in Section]?}, and the algorithm 
clearly captures the aftermath of this merger as only the purple and 
red tracer particles remain at z = 3. The filaments are highly per- 
turbed within the inner virial region and trace complex trajectories, 
no longer displaying the ordered gas inflow seen at higher redshift 
in the left panel. 



3.3 Locating the satellites 

There are two types of satellites in the simulation: ones with a dark 
matter halo and ones whose dark halo component has been stripped 
by tidal forces to the extent that it does not contain enough particles 
to be detected by the halo finder. One property all satellites share in 
common, however, is that they are substructures of the main host. 
As shown in the left panel of Fig. [3] luminous satellites in the simu- 
lation are often found to stream along filaments, eventually merging 
onto the central disk, perhaps transporting large amounts of angular 
momentum in the process. One is therefore faced with the follow- 
ing dilemma: should the angular momentum of gas gravitationally 
bound to a satellite, which itself is drifting along a filament, be 
solely attributed to the satellite? The answer to this question in this 
study is chosen to be 'no': only gaseous material within the virial 
region of a satellite that does not satisfy the nT criteria in that it lies 
above the upper filament density cut, is apportioned to the satel- 
lite phase. 'Satellite' is therefore made in reference to the central 
dense galaxy component which dominates the 'true' satellite angu- 
lar momentum signal. The advantages of this satellite definition are 
two-fold: 

• only small regions are masked from the density field, leaving 
the filament trajectories largely unperturbed, and; 

• the highest density material that is likely to dominate the satel- 
lite angular momentum budget is accounted for. 



The satellite identification procedure used in this study locates 
subhaloes of the main host that were once themselves field haloes 
earlier in their accretion history. Fig. [5] shows z projections of the 
gas density within the host virial region (orange circle) at z ~ 8 
(left panel) and z — 3 (right panel), and includes all the dark matter 
satellites (white circles) and stellar satellites (red circles) that have 
been found by applying the satellite identification technique, which 
is described in detail in Appendix[B| Each of these circles is centred 
on the centre of the galaxy it represents, as identified by the halo 
and stellar clump-finder. The filled black circles correspond to the 
constituent stellar particles of each stellar satellite, which can vio- 
late the virial accuracy condition imposed by the MSM halo-finder 
and can hence lie beyond their host's virial region, as illustrated 
for the elongated satellite in the right panel. Clearly the apparent 
luminous satellites in Fig.[3]are captured by the satellite-finding al- 
gorithm at the two extrema of the range in epochs examined in this 
study. 



4 ANGULAR MOMENTUM MEASUREMENT 
TECHNIQUES 

4.1 Calculating angular momentum in the general case 

Unless stated otherwise, the total angular momentum J has been 
computed by using the following expression: 

J = - r c ) x (u, - v c ), (4) 

i 

where i refers to either a stellar particle, a dark matter particle or 
a gas cell of mass rrii at location n with velocity vt, and r c and 
velocity v c refer to the centre of rotation within the central galaxy. 
There are two suitable estimators for these latter quantities: 1) the 
densest cell, which is likely to correspond to the disk centre, and; 
2) the centre of mass. Both of these definitions have been found 
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to yield very similar estimates of r c and v c in the NutCO run, but 
at high redshifts where the main host undergoes major mergers, it 
sometimes happens that the centre of mass is located in a region rel- 
atively devoid of mass, which does not correspond to the centre of 
the disk. We therefore adopt the densest cell as the best estimate of 
r c , which has been computed by simply summing the gas Q g , stel- 
lar Q-,, and dark matter 5dm grids within 0.1 r v h- of the host centre 
and locating the member cell of their sum C/o.i with the most mass. 
The velocity of the centre of disk rotation v c has been approxi- 
mated by the centre of mass velocity v conl of Go.i- Since there are 
a large C§> 10 6 ) number of cells and particles that contribute to 
Vcom at all times, noise is low and this velocity should thus be very 
close to the true bulk velocity of the disk. 

4.2 Obtaining the equation of the gas disk's plane of rotation 

It is assumed that the rotation of the gas disk is confined to a single 
plane centred on r c at all times, and an estimate for the equation of 
this plane has been made by using equation |4]( to compute the disk 
angular momentum J 4. All of the gas within C/0.1 associated with 
the luminous satellites (identified using the satellite-finding tech- 
nique) and the filaments (identified using the nT criteria) has been 
flagged, and the remaining gas in this region is assigned to the disk. 
Fig.|4]shows a density image of the gas disk at z = 3 in the NutCO 
run, and overplots the plane found by adopting the above disk def- 
inition. The excellent agreement between the computed disk plane 
and its apparent orientation in Fig.|4]provides confirmation that the 
disk definition and estimates of r c and v c in equation |4| are reli- 
able. 




Figure 4. A VAPOR visualization in density of the gas disk refined to 
I = 18 ( Ax ~ 12 pc) at z = 3, The black slice shows the plane of 
disk rotation found by computing the disk angular momentum within 
Go.l (equation]?). The scaled sidelength of the box in this image is approx- 
imately 6 kpc in physical coordinates, and the set of orthogonal axes and 
the units of the colour bar are the same as in Fig.[T] 

is populated by selecting all of the gas cells with a temperature 
above the imposed upper temperature limit of the cold phase (i.e. 
T > 2 x 10 5 K). Hot gas should not be regarded as a separate com- 
ponent, however, because these cells are included in the satellite 
and disk signals by construction (see Section |4~2} . 



4.3 Quantifying angular momentum transport to the disk 

The material channelled inward along filaments / in the simulation 
gets compressed once it settles into the disk d, and star formation 
ensues when the gas is sufficiently dense. Mass and angular mo- 
mentum transported by these flows is hence assumed to be shared 
amongst the central baryons b, as shown below by the arrows. The 
central baryons refer to both the gas in the disk and the stars that 
form from the disk's fragmentation: 



m/ - 
Jf 



m d + m* 
J d + J*. 



(5) 
(6) 



At a given time output, a list of identifiers of every particle that has 
ever been resolved as a member constituent of a satellite, as de- 
tected by the satellite-finding algorithm over the Milky Way host's 
virial region, is updated. Hence J* has been computed at each 
epoch by using equation (|4| and including each particle within 
0.1r v i r that does not belong to this list. These non-satellite stars 
correspond to bulge and disk stars, but the former component has a 
weak angular momentum modulus (see e.g. |van den Bosch et al.| 
|2002[ > and so the stellar disk signal dominates J*. Hence, even 
though m,b and Jb refer to the total (stars plus gas) mass and an- 
gular momentum of the 'baryonic galaxy' in what follows, we will 
often abbreviate this to 'the disk'. By excluding contributions from 
satellite stars, equations l|5j and j6]l focus on in situ disk formation. 

We now compare Jb at each time output with the accumu- 
lated angular momentum brought in by the satellite and cold fila- 
ment gas within 0.1r v i r , since these components are thought to fuel 
the disk's angular momentum budget. We also estimate the fraction 
of Jb acquired via a hot mode of gas accretion, as this enables di- 
rect comparisons with the cold mode signals. This hot gas phase 



4.3.1 The hot and cold phase contributions to mb and Jb 

Estimates of the angular momentum accumulated onto the disk 
via the accretion of hot and cold gas have been made by time- 
integrating the instantaneous flow rates of these components across 
a thin spherical shell located on the periphery of the disk region. 
This is a suitable method to employ because the inflowing gas in 
the simulation settles onto the disk very rapidly after crossing the 
shell (the validity of this assumption is tested in the following sec- 
tion). The flow rates have been computed by using the technique 
described in |Powell et al.| (201 1) , and for simplicity, reference is 
made to the filament gas phase in the descriptions below (but the 
same equations apply to the hot mode): 

dJ Ei PiVv.i [(ri - r c ) x (vi - v c )} (Axf 2 

"TT = — ^ , A x 3 47rr,„ (7) 

dt - E,(Ax)3 47rrm ' W 

The sum indexed by j is performed over all of the cells in the shell, 
irrespective of their gas phase, and r c and v c are the same as in 
equation |4|. The variables n,Vi and pi correspond respectively 
to the position, velocity and density of the i th cell satisfying the 
relevant imposed condition, which in this case is the requirement 
that i is a filament cell. Equations (|7} and {8} hence yield the net 
radial flux transported across the spherical shell towards (negative 
sign) or away from (positive sign) the disk, where the radial veloc- 
ity of each cell is given by: 



V r ,i = (Vi 



(9) 



It follows that the total angular momentum (and mass) transferred 
towards the disk by each filament over the time interval At = t — U 
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AJ flux (t) 




(10) 



(11) 



for the swept up gas, respectively. The ratio of these terms in equa- 
tion 1 12 1 satisfies AMfl ux /AAf s t a t > 10 for every time output in 



where the latter step approximates the integral by a simple numer- 
ical trapezium integration, i.e. a linear interpolation of the unre- 
solved flow rates between time outputs. The signal found using the 
'spherical flux technique' of equation \\\\ at t is then accumulated 
with the corresponding signals at earlier times. 

The summations indexed by i in equations {7} and {8} include 
all of the filament cells whose centres lie within a shell of inner, 
outer and mid radius r^ n = r m — Ar/2, r out = r m + Ar/2 and 
Tm = 0.1r v ir, respectively. The thickness of the shell Ar is a free 
parameter and has been fixed to 1% of r v i r (i.e. one tenth of the 
0.1r v i r radius) at each time output, as this choice strikes a good 
compromise between negligible Poisson error and the validity of 
the thin shell approximation (the shells typically included ~2 x 10 5 
cells at the highest redshifts around z ~ 9 and ~2 x 10 7 cells 
at the lower redshifts around z ~ 3, where the shell is physically 
larger). Similar flux measurements were found for shell thicknesses 
between ^0.1-2.5% of r v i r , which implies that there is some flex- 
ibility in the value assigned to Ar. A larger grid than Qo.i is there- 
fore required, in order to compensate for the half of the shell that 
spills beyond the 0.1r v i r region. The filaments are thus constructed 
over a 0.15r v i r region according to step 3 of the colour propaga- 
tion scheme described in Appendix [A] This particular grid length 
is chosen because the physical width of the cells within Q0.15 is 
still at the 12 pc limit for every time output, and so the extra grid 
length makes no difference to the accuracy of the mass and angular 
momentum computations. 

It is clear from the above that a decision has been made to 
represent the quantities r m and Ar as fractions of the host's virial 
radius, which grows with time in the simulation. The spherical flux 
method hence sweeps up gas by construction, and a simple trapez- 
ium integration of radial flux signals across the shell could miss 
this component if the filament's trajectory is mostly azimuthal. We 
therefore recalculate AMfl ux (and A Jfl ux ) to explicitly take the 
swept up component into account by: 

(i) fixing the sphere radius r(ti) over the given time interval 
At2i (= t% —tl) between time outputs t\ and ti\ 

(ii) applying equations {7]l and {8} to these identical spheres at 
the two time outputs, and performing the trapezium integration in 
equation \ W\ . and; 

(iii) combining the above spherical flux signals with the cor- 
responding instantaneous mass (AM sta t) and angular momentum 
(A Jstat) of all the filament cells measured at the later time output 
t<x within the small spherical shell spanning the volume between 
r(h) and rfa). 

The integrated quantities at time t hence correspond to: 

t 

M T(t) = ^[AM flux (^) + AM s tat(ii)] d 2 ) 



the simulation, indicating that the amount of gas swept up over a 
given time interval is much smaller than the radial mass flux (which 
is somewhat expected if the output sampling is fine enough). The 
initial time output in equations \12\ and {13} is denoted by ti and 
the sum is performed over all discrete time outputs i, ^ t\ ^ t, 
yielding projections along the disk rotation axis: 



J P {t) = J T (t) ■ J b (t). 



4.3.2 The satellite contributions to m b and Jt 



(14) 



If equation {8} were applied to the satellites, it could yield zero 
mass accretion events as satellite mergers are discrete burst events 
that can easily pass undetected through the shell over a given time 
interval. The following simple prescription, which seeks to pinpoint 
mergers with the disk, has been used to measure the instantaneous 
mass and angular momentum of gas belonging to satellites before 
it is accreted onto the central galaxy (satellite stars are ignored in 
this study as they are excluded from the baryonic disk definition 
in equation [5}. The algorithm begins by finding all the luminous 
satellites (according to the method described in Appendix [Bj that 
infringe the 0.1r v i r and r v i r regions at t and t + At respectively. 
Each satellite identification number at t is then compared with the 
list of satellite identifiers at t + At. Without a match, a given satel- 
lite has either fully merged with the disk or has been stripped of 
mass upon crossing the central region and has consequently passed 
below the minimum resolution limit of 100 particles imposed by 
the clump-finder. The merger has hence already taken place or is 
due to take place very shortly, and so the satellite signals at t are 
measured. If, however, the satellite is resolved at t + At (it may 
have moved closer to the disk or outside the central region during 
the merger process), then it has evidently not yet been accreted, and 
so no signal is recorded. The gas within the satellite virial region is 
hence included in the cumulative signal at t provided the satellite 
is: 

(i) encroaching the 0.1r v i r region at t, and; 

(ii) unresolved at t + At. 

The satellite gas mass M„ and angular momentum J s deposited 
within 0.1r v i r at t are subsequently computed using direct summa- 
tion and equation Q respectively, and are coupled with the totals 
from the previous epochs in an analogous fashion to equations l |12[ > 
and {13}: 

i 

AfJ«0.1r vi r(t),t) = ^M s «0.1r V i r (t;-),4) (15) 

*i 
t 

J^«0.1r vir (t),i) = ^Js« O.lTVir^),^), (16) 

*i 

yielding a projected component along the disk rotation axis: 
J"(t) = JtW- Jb(t). (17) 

The explicit time dependences of J p , Jt and J b in equations < | 1 4[ > 
and {T7} are henceforth dropped for brevity. 



Jt (t) = J2 1 A Jflux & ) + A Jstat & )] ' 



(13) 



where the first and second terms refer to the signals measured using 
the spherical flux method and the 'static shell' approximation used 



4.3.3 Assumptions 

By comparing at each epoch the projected angular momentum of 
gas measured using equations j I4[ i and {17} with the disk signal 
found using equation (16), two implicit assumptions are being made: 
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• the angular momentum of the material passing through the 
spherical shell at 0.1r v i r as measured at t reflects the amount de- 
posited onto the disk at the time of contact t + St (i.e. St is assumed 
to be smaller than the time interval between two consecutive out- 
puts of the simulation At), and; 

• the time integrated filament and hot gas flux signals are well 
approximated by a trapezium integration. 

The first condition draws attention to possible time lags between 
material crossing the shell and accreting onto the disk, and its im- 
portance can be estimated by comparing the free-fall time ts at 
0.1r v i r with the separation between consecutive time outputs in the 
simulation At (for a defense of the view that filaments are indeed 
well approximated by radial free-fall, see e.g. [Rosdahl &~ Blaizot 
|2012| l. In the NutCO run it was found that the range in ratio of 
these two timescales satisfied 0.1 < tg/At < 1, and so it is likely 
that the change in system configuration during the interval does not 
yield a significantly different amount of angular momentum trans- 
ferred to the disk compared with that measured at 0.1r v j r , thereby 
justifying the first assumption. The second assumption is consid- 
ered in the following section. 



5 RESULTS 

The techniques described in the previous section are now imple- 
mented in an attempt to examine whether: (i) filaments dominate 
Jb at high redshift, and; (ii) there are any special episodes of disk 
angular momentum acquisition. 

Fig.|5]plots the redshift evolution of the integrated net inflow 
of mass towards the disk in the upper panel, and the projection of 
the integrated angular momentum along the disk's axis of rotation 
in the lower panel. All of the measurements have been performed at 
or within the 0.1r v i r boundary as described in Section]?] Baryonic 
disk signals correspond to the solid black lines, while the satellite 
gas (equations 1 1 5 1 and 1 1 7} and hot gas (equations 1 1 2| and[T4]> sig- 
nals are given by the green and cyan lines respectively, with the 
filled circles in the lower panel indicating negative cumulative an- 
gular momentum projections. Each filament is distinguished by its 
unique colour, and the dashed dark blue line represents their sum at 
each epoch (equations|12|and[T4]>. The solid dark blue line shows 
the signals of all of the filament cells that have been found by ap- 
plying the nT criteria at each redshift, and serves as a consistency 
check for the individual filament computations. The dashed and 
solid dark blue lines are indistinguishable from one another in both 
panels, which implies that the vast majority of the filaments cells 
have been accounted for, and reinforces the high level of accuracy 
of the tracer colouring technique. Finally, it is expected that the 
total gas angular momentum and mass contributions from the fila- 
ment and satellite gas phases are equal to the baryonic disk signals 
at each epoch. The dashed black lines in Fig.[5]test this hypothesis 
and should be compared with the corresponding solid black lines. 

The disk, filament, hot gas and satellite signals at z' in ~ 8.5 
are subtracted from their counterparts at z < z' ln and are first 
recorded at z ~ 8 in Fig. [5] for reasons associated with the tracer 
colouring algorithm (which could not be applied until z- m ~ 10, 
the earliest time output that displayed a smooth, distinct filament 
configuration in the NutCO run), and the subsequent trapezium in- 
tegration of equations jT} and |8](. This approach of shifting the 
origin of the entire system forward to z' ln and modelling the 0.1r v i r 
sphere as being empty at z' in is justified given: a) the desire to main- 
tain consistency with the starting point of the tracer colouring al- 
gorithm, and; b) the negligible amount of time elapsed between z\ 



(the very first output from the simulation) and z[ n as compared to 
the time elapsed between z\ and the final output at z = 3. 

5.1 The mass growth of the central disk component 

It appears that the technical details surrounding the choice of start- 
ing point from the previous section are somewhat negligible be- 
cause the upper panel of Fig.[5]indicates that the baryonic disk ex- 
periences rapid growth at high redshift, where its mass increases by 
a factor of ~8 between z ~ 8 and z ~ 6. This growth is slower 
at later times, with a factor of ~2 increase in disk mass between 
z ~ 5 and z = 3. The periods of upturn in the green curves cor- 
respond to discrete luminous satellite merger events, which deposit 
~1% of the disk's mass in the form of gas at early times (z ~ 7) 
and ~10% at late times (z — 3). Satellites do not therefore ap- 
preciably affect the mass budget of the disk in the redshift range 
3 < z < 8. The amount of hot gas accreted onto the disk increases 
at a very stable rate owing to the long cooling times associated with 
a metal-poor intergalactic medium in the NutCO run, but by z — 3, 
the hot phase only accounts for about one tenth of the disk's accu- 
mulated mass. This phase thus plays a similarly subdominant role 
to the luminous satellites. Evidently the mass budget of the disk is 
largely controlled by the behaviour of filaments, and Fig. [5] sug- 
gests that it is the purple filament corresponding to region 1 in the 
left panel of Fig.[T]that carries the most mass to the disk, typically 
transporting at least twice as much mass at a given epoch compared 
with the yellow and red filaments, whose contributions flow in ap- 
proximately equal measure until z ~ 5.5 whereupon the purple and 
yellow filaments merge, with the former surviving (Fig.|2j. The dis- 
appearance of the yellow filament is evident from the plateau in its 
accumulated mass evolution for z < 5.5. 

By comparing the dashed black line with the solid black line 
in the upper panel of Fig. [5] it can be seen that the spherical flux 
method (which dominates the signal in equation [T2), designed to 
predict the net radial inflow of gas accreted onto the disk, conserves 
mass to high precision across the disk's evolution. This is remark- 
able given that the predicted flux rates across the shell could, in 
principle, vary quite dramatically between time intervals and hence 
not be well approximated by a trapezium integration. For exam- 
ple, it is possible that fast moving clumps of material that are ac- 
creted onto the disk pass through the thin shell undetected, which 
might lead one into thinking that the flux method has a system- 
atic tendency to underestimate the accumulated mass signal. The 
forepart of this argument appears to hold for the very first recorded 
mass signal at z ~ 8, where the inflow of gas is just under a fac- 
tor of 2 below the disk signal in place at that epoch. Clearly this 
is only a phenomenon that occurs at high redshift, as it quickly be- 
comes unimportant when the disk experiences rapid smooth growth 
via cold gas accretion over the subsequent time interval, which is 
captured by the flux method. Even so, for completeness, it is in- 
teresting to speculate on this high redshift discrepancy. One might 
propose that by increasing the thickness of the spherical shell from 
Ar — 1% of r v i r to a larger value, there is a higher probability 
of capturing this apparent 'missing' filament mass. Yet it has been 
found that even when Ar ~ 10% of r v i r is used at the first few time 
outputs, there is a negligible difference in the filament (and hence 
filament plus satellite gas) signals. The insensitivity of the flux rate 
to these changes at high redshift symbolizes the robustness of the 
spherical flux method and implies that the offset is not caused by an 
underestimation of the amount of mass transported radially inward 
by the filaments. Given the large 100 particle threshold imposed by 
the clump-finder and the technique adopted by the satellite-finding 
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Figure 5. The relative importance of hot gas (T > 2 X 10 5 K), cold filament gas and satellite gas in establishing the mass and angular momentum budget 
of the central disk at high redshift. Afj (upper panel) and Jt (lower panel) have b een accumulated with their corresponding signals from all prior epochs. 
The projection of along the unit disk direction J;, then yields J p (z) in equation jl4| for the filament and hot gas phases. A similar computation has been 
performed for the satellites, except that equations |T5j and \16\ have been used in place of equations \12\ and |T3j, with the filled circles corresponding to 
negative J p (z) contributions. The dashed dark blue lines sum over the individual filament signals at each epoch, and cannot be distinguished from the solid 
dark blue lines, which show the accumulated signals for all of the filament gas found by applying the nT cuts. The mass and angular momentum of the disk 
baryons defined by equations {5} and {6} have been computed at every time output and are given by the solid black lines. 



algorithm in Section [4.3. 2| of identifying satellites within 0.1r v i r , 
there does appear to be some scope for suggesting that the luminous 
satellite component is underestimated at this epoch: less massive, 
darting satellites could be excluded from the light green curves in 
Fig. [5] After all, the satellite contributions may appear somewhat 



artificial in the sense that it is unlikely that there are no satellite 
merger events beyond z ~ 7. We therefore surmise that this issue 
could be resolved by including the contribution of satellite galaxies 
coming from further away than 0.1r v i r at high redshift and pos- 
sibly reducing the particle threshold of the clump-finder, although 
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Figure 6. Top panel: A Mollweide projection of the trajectories followed by the cumulative angular momentum vectors of the disk and filaments as a function 
of redshift, colour coded as in Fig. [5] The filled circles indicate the starting point (z = 8) and the thick parts of the curves mark the redshift interval 
4-5 5; 2 < 5.5 where the cumulative angular momenta of the filaments undergo a dip in Fig. [5] Bottom panels: RGB images showing the z-projected 
temperature (red channel), density (green channel), and radial velocity (blue channel) of the gas enclosed within half the virial radius of the host halo during 
the redshift interval 4.5 < z < 5.5. Redshifts and scales are indicated at the bottom of each panel. Note the extremely disordered gas density field at z = 5.06 
close to the (smaller) disk in the middle bottom panel, as compared to the previous (bottom right panel) and later (bottom left panel) epochs. This disorderliness 
is correlated with both the minimum of the dip in filament angular momentum shown in Fig.[5]and the dramatic change in the direction of the angular momenta 
of all the filaments visible in the top panel above. 



this latter option would impact the reliability of our satellite mass 
and angular momentum measurements. In light of these difficul- 
ties, and the fact that this single interval of seeming lack of mass 
conservation constitutes the least relevant high redshift interval in 
Fig. [5] further analysis probing this discrepancy is deemed unnec- 
essary. Indeed, such an underestimate will only have a minor effect 
for z < 7.5 given the rapid rate of disk growth at early times. 



5.2 Which components dominate the disk's angular 
momentum budget? 

The general trends of satellites and hot gas contributing a weak 
mass signal are also imprinted in their angular momentum evolu- 
tion, shown in the lower panel of Fig. B] At z = 3, the hot gas 
signal (like the satellite signal) accounts for roughly one tenth of 
the disk's accumulated angular momentum, and with the exception 
of its evolution between 4.5 < z < 5.5 (discussed in more de- 
tail below), increases at a steady, uniform rate towards late times. 
In terms of angular momentum modulus, it therefore appears that 



the relative contributions of the satellite and hot gas components 
to the disk angular momentum budget generally trace their rela- 
tive disk mass contributions. In terms of angular momentum direc- 
tion, the satellite material which ends up deposited onto the central 
disk has a tendency to be anti-aligned with respect to the disk's an- 
gular momentum direction at high redshift, but because the signal 
strength is relatively weak, this mismatch only has a negligible ef- 
fect. Certainly at redshifts below z ~ 4.5, the mass that is headed 
towards the disk in the satellite (and hot gas) phase appears to be 
more closely aligned with the disk's direction. 

When integrated across all three components (solid and 
dashed dark blue lines), the filament angular momentum signals 
show a similarly strong correlation with their corresponding mass 
trajectories. Clearly for the redshift ranges 3 ^ z < 4 and 
5-5 < z < 7.5, the disk's budget is driven by the transport of angu- 
lar momentum from filaments. When decomposed into its individ- 
ual separate trajectories, however, the cold phase shows noticeably 
different behaviour: (i) the yellow and red filaments transport un- 
equal amounts of angular momentum along the disk's direction for 



12 H. Tillson et al. 



z > 4.5, in contrast with their relative mass contributions, and; (ii) 
the purple and yellow filaments are equally dominant for z > 5.5, 
unlike in the mass plot where the purple filament clearly dominates 
at all epochs. It should also be noted that the absence of a correla- 
tion between negative satellite signals and negative filament signals 
implies that satellites are not necessarily sensitive to the dynamics 
of the filament flow, despite often streaming along these flows (e.g. 
Fig# 

Arguably the most interesting feature of Fig.[5]is the lack of re- 
sponse the disk has to the sudden dip in the projected component of 
the purple filament's cumulative angular momentum during the fil- 
ament merger phase between 4.5 < z < 5.5. It was found that the 
change in the direction of the disk's axis of rotation during succes- 
sive intervals across this period was negligible (of the order a few 
degrees, as can be seen in the top panel of Fig.|6j, and so the sharp 
drop-off implies that the purple filament transports angular mo- 
mentum that is anti-aligned with the disk during the early phase of 
the merger (equation 1 14}. The red filament's disk-projected angu- 
lar momentum contribution also abruptly changes around z ~ 5.5, 
which when combined with the purple filament's, yields the flat- 
tened integrated dark blue curve that reflects the net small amounts 
of angular momentum being transported inwards by the filaments 
during this merger period. Looking at the bottom right and middle 
panels of Fig. [6] it becomes clear that this disruption is caused by 
the ongoing merger of the central disk with a quite massive satel- 
lite galaxy. This process leads to a complete reorganisation of the 
gas density field in the central region where the filaments connect 
to the disk, as shown in the bottom middle panel of Fig. [6] Whilst 
this does not affect the mass flow, as the gas must eventually pile 
up onto the disk in the absence of feedback and shock heating, it 
has dramatic consequences for the angular momentum of the gas 
embedded in the filaments. 

One may therefore expect the new population of disk stars that 
form from this anti-aligned angular momentum filament material to 
slightly reduce the disk angular momentum growth rate that is ob- 
served in the lower panel of Fig. [5] However, a new 'equilibrium' 
configuration between filaments and disk rapidly arises, presum- 
ably through torquing (see bottom left panel of Fig.[6](, in which the 
purple filament's angular momentum contribution along the disk's 
direction is restored. Moreover, at z < 4, the purple filament be- 
comes the dominant filament contributor to the angular momentum 
budget by at least a factor of ~3. This prominence is partly fuelled 
by the large-scale merger between the purple and yellow filaments 
that begins at z ~ 5.5 and induces a change in the mass, velocity 
and trajectory direction of the gas embedded within the surviving 
filament. Note that by z = 3, which corresponds to the end of the 
simulation, the angular momentum directions of the gas present in 
the red and purple filaments are fairly well aligned with one an- 
other but not with the disk (see the top panel of Fig. [SJ, which we 
interpret as small-scale torquing of large-scale driven quantities. 

Meanwhile, the angular momentum that our spherical flux 
method measures from filaments at 0.1r v i r differs from the fila- 
ment angular momentum deposited onto the disk during the satel- 
lite merger period because the inner dynamics within 0.1r v i r are 
severely perturbed (note that the resultant complex inner gas tra- 
jectories appear to be imprinted in the hot gas signal too, which be- 
comes anti-aligned in Fig.[5]around z ~ 5). It is likely that isother- 
mal shocks in the disk boundary region are responsible for filtering 
aligned angular momentum to the disk during this phase, driving its 
angular momentum growth, because the baryonic disk angular mo- 
mentum is dominated by its stellar component (| J*|/| J gas | ~ 2.6 
at z ~ 5 for example), which seems to stabilise it against strong lo- 



cal perturbations. Should this effect apply to a large fraction of disk 
galaxies, one could conclude that contrary to dark matter haloes 
(e.g. |Bett & Frenk|2012) , minor mergers of galaxies do not seem 
capable of significantly altering the direction of the angular mo- 
mentum of the central disk. 



6 DISCUSSION 

Perhaps the most defining feature of disks in galaxies with M* < 
3 x 10 10 Mq is the amount of angular momentum they possess, as 
it is this property that establishes many of their fundamental scaling 
relations such as the Tully-Fisher relation (e.g. Governato et al. 
|2007[ > and those linking velocity, metallicity, surface brightness and 
stellar mass jDekel & Woopfe] |Dekel & Birnboim|[200o) . The 
Tully-Fisher relation, in particular, is an important constraint that 
many semi-analytic models of galaxy formation strive to satisfy 
( |Hatton et al.pM] |Cattaneo et al.|2006| |Croton et al.p006l [De] 
|Lucia et al.|2006| l. Understanding how the disk acquires its angular 
momentum is hence of central importance, and has been the subject 
of this paper. 

6.1 The dominance of filaments at high redshift 

The evidence favouring the growth of disks in low mass systems via 
a mode of cold gas accretion appears to be mounting. Keres et al. 
( 2005 1 performed SPH simulations of several hundred galaxies be- 
tween ^ z ^ 3 using a comoving gravitational softening scale of 
5 kpc, and found evidence of a clear shift in the total fraction 
of gas accreted with temperatures T < 2.5 x 10 5 K around a galaxy 
stellar mass of ~ 2.5 x 10 10 Mq, with at most a factor of ~1.6 
deviation in M* across this redshift range. |Ocvirk et al.] < [2008) > con- 
ducted a similar statistical analysis on haloes from the HORIZON- 
MARENOSTRUM simulation with 10 10 < M H /M < 10 13 at 
higher redshifts between 2 sC z 5, in an attempt to measure 
the temperature and density of the gas deposited onto the central 
galaxies of these systems, which were simulated using the RAM- 
SES code at a physical resolution scale of ~ Ih^ 1 kpc. Their re- 
sults demonstrated that the fraction of gas accreted at temperatures 
T < 2.5 x 10 5 K sharply increases for M H < 4 x 10 11 M , 
in concordance with the estimates provided by Birnb oim & Dekel] 
(2003 ) and Keres et al. (2005 1. By resolving and examining the var- 
ious components of the filamentary gas phase across a wide range 
in redshift, the results from Fig.|5]in this paper provide quantitive 
support for the cold mode paradigm of gas accretion onto a Milky 
Way like disk, advocated by the above studies. They also show that 
a single filament can be responsible for driving the mass budget of 
the baryonic disk, although this does not necessarily map to a dom- 
inance of the disk's angular momentum budget, as two filaments 
appear to closely share priority at the higher redshifts (z > 5.5), 
despite a clear difference in their mass contributions. The flow of 
cold gas onto the central disk is not necessarily a smooth, continu- 
ous process either. While satellites do not significantly perturb the 
filament trajectories in general (except in the region close to the 
central disk, and even then, not very often and only over short pe- 
riods of time), large-scale motions of filaments can lead to mergers 
between these cold streams jPichon et al.|20lT| and this process 
is able to change the orientation of angular momentum advected 
along these flows, probably even if gas channelled along each of 
the filament mergers in question is co-planar with the disk's rota- 
tion before the merger. 

A quantitative result of this kind has remained elusive to many 
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previous studies mostly due to the parsec-scale physical resolution 
required to resolve the satellites, filaments and disk components 
within the inner tenth of the virial region at z > 3. Some au- 
thors have recently started to probe this resolution barrier, however. 
Kimm et al. (201 1 1 measured the angular momentum of gas in ra- 
dial bins as a function of redshift, for the same resimulated NUT 
Milky Way like halo analyzed in this paper. Upon stacking the ra- 
dial profiles in two separate redshift regimes (z ^ 3 and z > 3), 
they found evidence for a sudden loss in the amount of specific an- 
gular momentum transported by gas in the 0.1r v ; r region, but with 
a physical resolution scale of only ~50 pc, were unable to speculate 
on the cause of this loss. Although the analysis in this paper does 
not address this issue per se, it does show that the filaments supply 
enough angular momentum to the boundary of the disk region to 
be able to account for the angular momentum locked-up in the disk 
across most epochs at z > 3, a result that has been found by using 
a filament tracer colouring technique. 

Several studies have also demonstrated that the angular mo- 
mentum vectors of dark matter and gas are not necessarily aligned 
within the virial region of galaxies at both low (z — 0) and high 
(2 sC 3) redshift ( TBett et al.|20Tol |Roskar et al.|2010| >. |Danovich| 
|et al.| |2012l reported a weak correlation between the angular mo- 
mentum direction of gas in the inner (~0.1r v i r ) and outer (~lr v i r ) 
galaxy regions at z = 2.5, and a large body of evidence has now 
been presented that confirms the existence of mismatches in gas an- 
gular momentum directions over scales of a few disk scalelengths, 
giving rise to the population of 'warped' disks iGarc ia-Ruiz et al] 
|2002| [Shen & Se llwood 20061. Despite being observed as a low 
redshift phenomena, the recent simulation study by Rosk ar et al.| 
( |2010} has argued that warps may exist around z ~ 2-3. In Fig.[5]of 
this paper, the amount of angular momentum transported along the 
disk's direction by all the filaments is comparable to that locked-up 
in the disk, which would not be the case if large-scale misaligned 
signals reported above were preserved at the disk boundary. If the 
claim that the directions of the gas angular momentum vectors vary 
with distance from the central galaxy holds across the galaxy's life- 
time, some physical mechanism must be responsible for aligning 
the infalling angular momentum at the disk's edge. For example, 
Roskar et al. (20101 have argued that freshly accreted gas at the 
virial radius is strongly torqued by the hot halo gas component, and 
cited this as a possible cause of the warps between inner and outer 
disk structure. Whilst hot gas torques are unlikely to be dominant in 
our case, it is nevertheless surmised that shocks and/or torques re- 
move the components of the purple filament's velocity that are mis- 
aligned with respect to the disk's plane of rotation, thereby explain- 
ing why this filament is able to maintain a relatively stable transfer 
of aligned angular momentum to the disk between 3 < z < 8. 

Despite the evidence from both Fig.|5]and the aforementioned 
simulation studies favouring the cold gas paradigm, not all authors 
are convinced that this phase is quite so dominant in growing the 
disk components of low mass galaxies. Murante et al. (2012| mon- 
itored the fractional accretion rate onto two Milky Way like haloes 
and found that ~50% of the accreted gas onto the central galaxy 
between 3 ^ z ^ 6 was in a warm phase, with a temperature 
in the range 2.5 x 10 5 < T/K < 10 6 . They attributed this ap- 
parent discrepancy to a supernovae feedback prescription that only 
modelled thermal heating, as opposed to thermal and kinetic heat- 
ing. The effects of supernovae feedback on filament structure were 
also examined by |van de Voort et al.| ( [20lT| l, who suggested that 
it reduces the growth rate of low mass central galaxies residing in 
haloes with mass A/h < 10 12 Mq, implying that the filaments 
streaming cold gas at large inflow rates in Fig. [5] do not survive 



when supernovae feedback is included. |Powell et al.| ((201 1 ), how- 
ever, argued that this result is probably an artefact of supernovae 
feedback implementation, because with individual Sedov-Taylor 
blasts from the ultra-high resolution (Ai res ~ 0.5 pc) NUT feed- 
back run simulation, the net mass inflow rates of gas in the fila- 
ment phase were found to be an order of magnitude larger than 
the supernovae-driven mass outflow rates. The |Powell et al.|{201 1| > 
study hence implies that whilst the spatial distribution of filaments 
local to the disk is likely to be perturbed around sites of intensive 
stellar explosions, the amount of filament gas flowing towards the 
disk should largely be unaltered. Therefore, the amount of angular 
momentum aligned with the disk that is transported by filaments 
is probably unchanged too, given the general correlation between 
the filament mass and disk-projected angular momentum signals in 
Fig. [5] It is therefore expected that the difference between the su- 
pernovae feedback version of Fig. [5] and its counterpart presented 
in this paper is not significant, but we leave this comparison as a 
future exercise. 



6.2 Understanding the impact of mergers 

Fig. [5] shows that luminous satellite mergers constitute only a mi- 
nor fraction of the disk's mass and angular momentum modulus 
between 3 < z < 8. This result does not necessarily imply, how- 
ever, that satellites as a general population have little effect on the 
disk's evolution. |Bett & Frenk| ( [20T2 1 analyzed present day Milky 



Way haloes of mass M H ~ lO^-lO 12 5 ft" 1 M from one of the 
Millennium Simulation runs and examined the importance of 'spin 
flips', which are defined as abrupt changes of more than 45° in 
the orientation of a component's angular momentum vector. They 
argued that over 90% of their detected host halo spin flips were 
caused by minor mergers, and demonstrated that the number of flip 
events increases in the inner halo where the central galaxy resides. 
They further speculated that these spin flips could destroy the host's 
stellar disk (which is included in the baryonic disk signal of Fig.|5]l, 
or torque it (see, for example, Ostriker & Binney 19891. In both of 
these scenarios, one would expect a change in the angular momen- 
tum direction of the baryonic disk component. The validity of this 
hypothesis remains an open question because the existence of a 
simple correlation between dark halo spin flips and disk spin flips 
is yet to be confirmed: several recent studies have in fact hinted 
at a lack of correlation between the two ([Scannapieco et al.|2009{ 
|Stinson et al.|2010|[Sales et al.|2011J. It will be in teresting to test 
in more detail whether this claim by Bett & Frenk ( 2012 1 holds for 
the NUT host halo in this study, which experiences multiple minor 
mergers across its accretion history (as indicated by Fig. [3}. How- 
ever, we suggest that it is highly unlikely since the top panel of 
Fig. [6] clearly shows that the angular momentum direction of the 
disk varies very little throughout the majority of its lifetime. 



6.3 High redshift contributions to the present day mass and 
angular momentum of the disk from the NutCO run 

It is also informative to estimate the fraction of the disk's mass and 
angular momentum at z — that is already in place at z = 3, 
as large high redshift contributions would further highlight the im- 
portance of understanding the primordial phase of Milky Way like 
galaxy growth. Crude estimates are hence provided in this section, 
but it should be noted that more rigorous measurements could be 
made by analyzing all of the lower resolution outputs from the 
NutCO run between ^ z ^ 3 and extending Fig.[5]to 2 = 0. 
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|Kimm et al.| ( |2011[ ) measured the total amount of baryonic 
mass (gas plus stars) within 0.1r v i r from the Ai m i n = 48 pc 
NutCO ran at z = and found that M b (z = 0) ~ 8 x 10 10 M . 
This measurement includes the central contributions from all of the 
gas phases (satellite, hot mode, cold mode and disk) and all of the 
stars (satellite, disk and bulge), and is a factor of ~6 greater than 
the baryonic disk signal at z = 3 from Fig. [5] which is Mb(z = 
3) ~ 1.3 X 10 10 Mq. Multiplying the modulus of the total bary- 
onic specific angular momentum contributions within the O.lrvir 
region at z = from the Kimm et aE| ( |201 1] > study by the relevant 
component masses yields an estimate of the total baryonic angu- 
lar momentum signal: J b {z = 0) ~ 4.4 x 10 13 Mq kpckms" 1 . 
This is a factor of -~30 higher than the modulus of the bary- 
onic disk angular momentum signal at z = 3 shown in Fig. [5] 
J b (z = 3) ~ 1.6 x 10 12 Mq kpckms" 1 . Note that the mass 
(and most probably the angular momentum) ratios are upper limit 
estimates because the Kimm et al. (20111 measurements include 
the satellite star and non-disk gas contributions, whereas Fig.[5]just 
shows the baryonic disk signals defined according to equations |5j 
and {6). 

In summary, at least ~16% (~4%) of the NUT disk's final 
mass (angular momentum) is already in place by z — 3, suggest- 
ing that high redshift epochs studied in this paper represent a non- 
negligible period of the disk's accretion history. 



7 CONCLUSIONS 

We have analyzed an Adaptive Mesh Refinement cosmological res- 
imulation with 12 pc resolution in order to address whether the an- 
gular momentum acquired by the baryonic disk of a Milky Way like 
galaxy at z ^ 3 is driven by filaments. The filaments have been 
identified by tagging Lagrangian tracer particles within a region 
of length 27Vj r from the disk's centre at each epoch. The method 
assigns one of three colours (one for each filament) to the tracer 
particles by using a nearest neighbours scheme, and subsequently 
locates the 'pure' cells within the spherical grid of radius 2r v j r that 
contain tracer particles of identical colour. Individual filament tra- 
jectories are then grown around these pure sites. A satellite finding 
algorithm has also been presented, which operates under the princi- 
ple that a given satellite substructure of the main host at time t has 
at some earlier time been a separate field galaxy of its own. With the 
filaments and satellites resolved, the aim of this paper has been to 
measure their contributions to the host's disk angular momentum 
budget at each epoch. For the filament and hot gas components, 
the mass and angular momentum transported to the central disk has 
been computed by recording both the net inward flux across an out- 
wardly propagating spherical shell on the disk's edge at 0.1r v i r , and 
the contribution from the material swept up by the shell's outward 
radial motion. Satellite signals have been measured by adopting 
an alternative technique that pinpoints mergers with the disk and 
records their angular momentum at the final pre-accretion stage. 
The results suggest that: 

• The cumulative effect of the cold filament gas phase domi- 
nates the mass and angular momentum budgets of the disk as a 
function of time, hence providing quantitative support for the cold 
gas paradigm of disk growth in low mass galaxies at high redshift. 

• For 5.5 < z < 7.5, the largest portion of the filament angu- 
lar momentum signal is transported by two filaments, which start 
to merge around z ~ 5.5. This is accompanied by the central disk 
merging with a satellite galaxy, which temporarily induces a dra- 
matic change in the orientation of the angular momentum trans- 



ported towards the disk, but the surviving filament remnant quickly 
realigns itself along the direction imposed by the large-scale flow 
and by z ~ 4 dominates the filament angular momentum budget 
by at least a factor of 3. 

• The luminous satellites account for at most one tenth of the 
disk mass and angular momentum modulus at any given epoch. 
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APPENDIX A: IDENTIFYING THE FILAMENTS USING 
TRACER PARTICLES 

Al The algorithm 

The tracer propagation method, first encountered in Section |3?2"] 
consists of three principal stages: 

(i) Colouring the tracer particles across time. At some initial 
early time ti (zi ~ 10), planar cuts are used to divide the filaments 
into distinct regions within a spherical grid Q2.0 of radius 2r v - 1T 
(black planes in the left panel of Fig. [TJ. Each region is assigned a 
unique colour and tracer particles are initially coloured according 
to the region in which they reside. Hence there is a single filament 
per region by construction, and we refer to the filaments in regions 
1, 2, and 3 of the left panel of Fig. [T] as the purple, yellow and red 
filaments respectively. 

The Lagrangian nature of the system poses two computational 
difficulties when assigning colour: 

• filaments mix in the vicinity of the disk (as can be seen in 
the right panel of Fig.[TJ, and; 

• filaments acquire a drift velocity on scales r > r v i r ( |Pichon| 
|etal.|2011) , which can lead to them switching between regions 
within r v i r as the host grows. 

The first issue is tackled by resolving the filament trajectories at 
the following time output £2, by which time the vast majority of 
the incorrectly coloured particles within 0.1r v i r at ti have accreted 
onto the disk during the At = ti — ti interval, and are therefore no 
longer associated with a filament. To circumvent the second issue, 
previously flagged tracer particles i at locations ri(t\) that have 
moved during the At interval are updated according to their new 
instantaneous locations Ti{t2) at the following time t%, provided 
that they coincide with filament cells at ti . All of the newly accreted 
tracer particles that are co-spatial with filaments within the (larger) 
host spherical region at t-z that have not been assigned a colour are 
then given a colour according to their host region. At subsequent 
time outputs tj (£3 onwards), the newly accreted tracer particles 
are instead assigned the majority colour of their TV nearest coloured 
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Figure Al. These VAPOR images illustrate the trajectories of each gas filament on So. 5 (r v ; r ~ 7.5 kpc) at z ~ 8 in the NutCO run with i = 17 (Ax ~ 
12 pc), as identified by the colour propagation technique. Each panel corresponds to the updated versions of the purple, yellow and red filaments associated 
with regions 1, 2 and 3 respectively in the left panel of Fig.[T] and follows the same density colour scaling as in previous figures. The black curved arrow in 
the top panel indicates the gas flow direction of the purple filament. 



tracer neighbours that sample the filaments at r i (t 3 ■) (N was at least 
100). Note that the above neighbouring scheme cannot be used for 
t2 because the constructed filament trajectories are not known until 
the end of the time output (see steps 2 and 3 below). 

(ii) Colouring the grid. Once all of the filament tracer particles 
within Q2.0 are coloured, the overall colours of the member tracer 
cells are determined (particles located in non-filament cells within 
r v ir are ignored, so that only filaments are resolved). The prescrip- 
tion used in this study enforces the condition that each coloured 
tracer cell is 'pure', consisting entirely of particles with identical 
colour. The advantage of this approach is that it avoids the use of 
ad hoc criteria to help decide the overall colour of a cell that con- 
tains a certain ratio of particles of different colours. 

(iii) Individual filament growth. The final step propagates the 
pure cell colours to neighbouring uncoloured filament cells within 
C/2.0- Filaments are simultaneously grown around the cells with 
pure colour and the procedure terminates once all neighbours of 
neighbours etc, are coloured and the only remaining neighbouring 
cells are non-filament cells or cells that already have a colour. 

A2 Uncertainties 

While the nT cuts are only approximate, |Powell et al.| l |2011) have 
argued that the associated uncertainties in the bulk physical proper- 



ties of the filaments (average mass, velocity etc) remain very small 
because: (i) these cuts successfully capture the dense core of the 
filaments (as confirmed by visual verification) around the disk in 
the NUT suite runs with Ax rcs = 0.5 pc, and; (ii) reducing the 
lower density threshold has little impact on these properties. Hence 
we also ignore the uncertainties associated with the nT criteria in 
what follows, and the constructed filaments are calibrated against 
their nT counterparts. 

The uncertainties introduced by the filament colour propaga- 
tion scheme can be of two types: 

(i) not attributing any colour to a valid filament cell, and; 

(ii) attributing the wrong colour to a filament cell. 

The relative importance of uncertainties of type (i) was assessed 
by recording the number of filament cells without a colour upon 
application of the propagation method, as a fraction of the total 
number of filament cells found by applying the nT criteria, for 
different grid sizes at z ^ 3 in the NutCO run. These uncertainties 
were measured to be below the percent level at all times. 

Type (ii) uncertainties are now examined. Since the aim of 
this paper is to quantify the amount of angular momentum trans- 
ferred to the disk by filaments, it is key that the trajectories of these 
cold flows in the disk's vicinity are reconstructed to a high level of 
precision. The filament reconstruction stage (step 3 of the colour 
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propagation algorithm) was hence repeated over the 0.15r v i r re- 
gion to ensure that the angular momentum calculations around the 
disk were performed at the resolution limit (this particular choice 
of grid radius is justified in Section [4. 3. 1} . Once constructed, it is 
natural to wonder whether the filament trajectories are representa- 
tive of their 'true' counterparts in the simulation. A rigorous quan- 
tification of the relative discrepancies has proven elusive, so we 
provide the reader with a powerful qualitative diagnostic instead: 
a computer assisted visual check where the filament density field 
identified by the tracer colouring algorithm is compared to calcu- 
lated flow lines. Fig. |A1| shows 3D images of the density of each 
filament from the right panel of Fig. [T] as found by implementing 
the colour propagation technique on C/o.5 with spacing Ax ~ 12 pc 
at z ~ 8. By zooming into the larger 0.5r v i r region at the resolu- 
tion limit, it is possible to test the algorithm's ability to follow each 
filament onto the disk. The viewing angle in each panel of Fig. |A1| 
has been rotated to demonstrate the nature of the individual trajec- 
tories, and varies from filament to filament. By comparison with 
the right panel of Fig.[TJ it is clear from Fig. | A 1 1 that the separate 
filament paths are correctly identified to the point of disk contact, 
and that the curvature about the disk region is particularly well- 
resolved, which is most evident for the yellow filament in region 
2. Note that the sudden upturn of the red filament's inner trajectory 
in Fig. [Tj after mixing in the central region is also captured by the 
colour propagation routine. Thus Fig. |Al| strongly suggests that the 
algorithm is able to successfully construct the trajectories of indi- 
vidual filaments even though, contrary to type (i) uncertainties, it is 
difficult to argue that it does so at the percent level accuracy. 



APPENDIX B: PINPOINTING THE SATELLITES 

The details of the satellite-finding algorithm introduced in Section 
|3.3| are provided in this Appendix. 

In this study, a substructure i that is classed as a satellite of the 
host at a given time output ti satisfies both of the following criteria: 

• it partially or fully infringes the host's virial sphere at ti, and; 

• it has previously been a halo. 

The former condition ensures that only substructures within the 
host's virial region are considered, and is hence satisfied should 
any region within satellite i's virial radius also be contained by the 
host's virial sphere. If it is, then i is flagged as a satellite candidate, 
otherwise it is skipped. Assuming i is a satellite candidate, the lat- 
ter condition involves computing its merger tree, and if at any time 
output along its main branch one of its main progenitors becomes 
a halo, then i is flagged as a satellite. 

In order to compute the main branch along i's merger tree, 
the progenitor that donates the most mass to i is identified: this 
progenitor, j, is called the 'main progenitor', and the link from i 
to j is referred to as the 'main branch' (following Tillso n~et al.| 
|2011) . The satellite-finding algorithm imposes the condition that in 
order for a subhalo to be a satellite at time t, it must have been a 
field halo at some point earlier in its merger history. Note that this 
method allows objects to share the same main progenitor and hence 
maximizes the number of possible satellites, resulting in the purest 
filamentary angular momentum signals once the central satellite 
galaxy regions have been masked according to the prescription in 
Sectionl33l 

This paper has been typeset from a TpX/ DTpX file prepared by the 
author. 



